Newtons Method 

Newton-Raphson method or only Newton's method is one of the more popular methods used for solving nonlinear algebraic equations of the form f(x) = 0. Nonlinear algebraic equations contain powers of variable(s) and/or transcendental functions. One of the major achievements in mathematics was the proof (by Niels Henrik Abel (1802-1829) in 1824 that polynomial equations of degree greater than 4 cannot be solved by means of an algebraic formula (that is, in terms of radicals). 

The basis for Newton's method is that the actual root is estimated and the zero of the tangent to the function at that point is determined. 

Image 

 

If a real root x[0] is to be assumed for the function f, then one may easily compute fx[1]If we draw a line tangent to the curve at fx[0] then the tangent line intersects the x-axis at a point which is expected to be closer to the actual root r than the assumed root From the figure we see that  

 

solving for x[1]gives 

 

If we repeat the procedure at fx[1]and let x[2]be the x-intercept of the second tangent line, we find 

 

Repetitive use of the procedure yields a sequence of approximations that we expect converges to the root  

 

Newton-Raphson Method 

Suppose f is differentiable and suppose r represents the unknown root of f(x) = 0.Let x[0] denote a number that is chosen arbitrarily as a first guess to then repetitive use or iteration of  

 

yields a sequence x[1], x[2], x[3] .. () of approximations that we expect converges to the root  

Newton-Raphson Method 

 

Note that Newton's method unfortunately does not always converge to a root if x[0]is not sufficiently close to the actual root. 

 

 

 

 

See also NewtonsMethod . 

 

The formula of a fourth-degree equation is quite complicated as Maple shows below. 

> restart:MathMaple:-ini():
 

> f:=x->x^4-x-1;
 

proc (x) options operator, arrow; `+`(`*`(`^`(x, 4)), `-`(x), `-`(1)) end proc
 

> solve(f(x)=0,{x}):
 

> evalf(%);   #evaluates to floating points
 

{x = 1.220744085}, {x = `+`(`-`(.2481260628), `*`(1.033982061, `*`(I)))}, {x = -.7244919590}, {x = `+`(`-`(.2481260628), `-`(`*`(1.033982061, `*`(I))))}
{x = 1.220744085}, {x = `+`(`-`(.2481260628), `*`(1.033982061, `*`(I)))}, {x = -.7244919590}, {x = `+`(`-`(.2481260628), `-`(`*`(1.033982061, `*`(I))))}
 

The roots of the equation above can of course be approximated numerically using Maple's   fsolve . 

> fsolve(f(x)=0,{x});   #all real values  
 

{x = -.7244919590}, {x = 1.220744085}
 

> fsolve(f(x)=0,{x},complex); #all values
 

{x = -.7244919590}, {x = `+`(`-`(.2481260628), `-`(`*`(1.033982061, `*`(I))))}, {x = `+`(`-`(.2481260628), `*`(1.033982061, `*`(I)))}, {x = 1.220744085}
{x = -.7244919590}, {x = `+`(`-`(.2481260628), `-`(`*`(1.033982061, `*`(I))))}, {x = `+`(`-`(.2481260628), `*`(1.033982061, `*`(I)))}, {x = 1.220744085}
 

 

The solution of  

> f:=x->x^5-3*x^3+x^2-23*x+19:
f(x)=0;
 

`+`(`*`(`^`(x, 5)), `-`(`*`(3, `*`(`^`(x, 3)))), `*`(`^`(x, 2)), `-`(`*`(23, `*`(x))), 19) = 0
 

can be found using fsolve . 

> fsolve(f(x)=0,{x},complex);
 

{x = -2.722493355}, {x = `+`(`-`(.1943723329), `-`(`*`(1.932038154, `*`(I))))}, {x = `+`(`-`(.1943723329), `*`(1.932038154, `*`(I)))}, {x = .8012614801}, {x = 2.309976541}
{x = -2.722493355}, {x = `+`(`-`(.1943723329), `-`(`*`(1.932038154, `*`(I))))}, {x = `+`(`-`(.1943723329), `*`(1.932038154, `*`(I)))}, {x = .8012614801}, {x = 2.309976541}
 

> fsolve(f(x)=0,{x}); #real solutions
 

{x = -2.722493355}, {x = .8012614801}, {x = 2.309976541}
 

> plot(f(x),x=-3..3,`f(x)`);
 

Plot_2d
 

 

The following shows an iteration process which does not converge. 

> f:=x->piecewise(x>0,sqrt(x),-sqrt(-x)):
'f(x)'=f(x);
 

f(x) = piecewise(`<`(0, x), `*`(`^`(x, `/`(1, 2))), `+`(`-`(`*`(`^`(`+`(`-`(x)), `/`(1, 2))))))
 

First using Newton  

> Newton(f,1,10,0.0001);
 

-1.000000000, 1.000000000, -1.000000000, 1.000000000, -1.000000000, 1.000000000, -1.000000000, 1.000000000, -1.000000000, 1.000000000
-1.000000000, 1.000000000, -1.000000000, 1.000000000, -1.000000000, 1.000000000, -1.000000000, 1.000000000, -1.000000000, 1.000000000
 

In this case Newton's method fails. The animation below indicates the reason. 

> NewtonPlot(f, 1, 10,0.0001,grey);
 

Plot_2d
 

 

Example 1 

a) The volume of a body is V(x) = `+`(`*`(3, `*`(`^`(x, 3))), `-`(`*`(16, `*`(`^`(x, 2)))), `*`(20, `*`(x)), `-`(40)). Find x by solving the equation V(x) = 60. 

b) Solve `+`(`*`(18, `*`(sin, `*`(x))), `-`(`*`(6, `*`(sqrt(`+`(x, 1)))))) = 0. 

Solution 

> restart:MathMaple:-ini():
 

a) 

> V:=x->3*x^3-16*x^2+20*x-40;
 

proc (x) options operator, arrow; `+`(`*`(3, `*`(`^`(x, 3))), `-`(`*`(16, `*`(`^`(x, 2)))), `*`(20, `*`(x)), `-`(40)) end proc
 

> plot(V(x)-60,x=2..6);
 

Plot_2d
 

We use the formula where f(x) = `+`(V(x), `-`(60.)) 

> f:=x->V(x)-60:
 

We start start with x(0) = 5 and use 6 iterations. 

> x(0),N:=5,6;
 

5, 6
 

> for n from 0 to N do
   x(n+1):=evalf(x(n)-f(x(n))/D(f)(x(n)));
od;
 

 

 

 

 

 

 

5.294117647
5.268981078
5.268784059
5.268784047
5.268784047
5.268784047
5.268784047
 

The approximate value is  

> x:='x':
 

b) 

> f:=x->18*sin(x) - 6*sqrt(x+1):
f(x)=0;
 

`+`(`*`(18, `*`(sin(x))), `-`(`*`(6, `*`(`^`(`+`(x, 1), `/`(1, 2)))))) = 0
 

> plot(f(x),x=0..15,tickmarks=[10,10]);
 

Plot_2d
 

The figure shows that there is 3 or 4 roots. 

Newton method gives 

x[0] = 1 

> Newton(f,1,6,0.00001);
 

[.1240018122, .3990677339, .4063891309, .4064008634, .4064008631]
[.1240018122, .3990677339, .4063891309, .4064008634, .4064008631]
 

The solution is  

x[0] = 2 

> Newton(f,2,6,0.00001);
 

[2.647863728, 2.479722108, 2.471515518, 2.471492230, 2.471492231]
[2.647863728, 2.479722108, 2.471515518, 2.471492230, 2.471492231]
 

The solution is  

x[0] = 7 

> Newton(f,7,6,0.00001);
 

[7.411269108, 7.581510032, 7.644488093, 7.657140144, 7.657702263, 7.657703379]
[7.411269108, 7.581510032, 7.644488093, 7.657140144, 7.657702263, 7.657703379]
 

The solution is  

x[0] = 8 

> Newton(f,8,6,0.00001);
 

[7.947070592, 7.937756334, 7.937447443, 7.937447100]
[7.947070592, 7.937756334, 7.937447443, 7.937447100]
 

The solution is  

Using fsolve , we get 

> x[1]=fsolve(f(x)=0,x=1);
 

x[1] = .4064008633
 

> x[2]=fsolve(f(x)=0,x=2);
 

x[2] = 2.471492230
 

> x[3]=fsolve(f(x)=0,x=7);
 

x[3] = 7.657703381
 

> x[4]=fsolve(f(x)=0,x=8);
 

x[4] = 7.937447102
 

To go from one approximation to the next approximation can be visualized by animation. 

> NewtonPlot(f,7,6,0.00001,yellow);   
 

Plot_2d
 

Click on the figure and start the animation. 

> NewtonPlot(f,8,6,0.00001,grey);
 

Plot_2d
 

Let us start with a larger initial value. 

> NewtonPlot(f,12,10,0.00001,grey);
 

Plot_2d
 

Click on the figure and select Next Frame repeatedly.